import argparse
import numpy as np

def process_file(file_path, num_lines):
    signal_means = []
    single_base_signal_lengths = []
    seq_lengths = []
    mapping_quality = []
    shift_dacs_to_pa = []
    scale_dacs_to_pa = []
    shift_pa_to_norm = []
    scale_pa_to_norm = []

    with open(file_path, 'r') as file:
        for i, line in enumerate(file):
            if num_lines > 0 and i >= num_lines:
                break
            
            columns = line.strip().split('\t')
            if len(columns) < 15:
                continue  # 跳过无效行
            
            # 计算信号的均值
            norm_signals_text = columns[4]  # 第5列 norm_signals_text
            signal_segments = [
                list(map(float, segment.split(',')))
                for segment in norm_signals_text.split(';')
                if segment.strip()  # 跳过空信号
            ]
            
            # 每行信号均值
            signal_mean = np.mean([value for segment in signal_segments for value in segment]) if signal_segments else 0
            signal_means.append(signal_mean)
            
            # 单个碱基信号长度
            single_base_signal_lengths.extend([len(segment) for segment in signal_segments])
            
            # 计算 seq 的长度
            seq = columns[3]  # 第4列 seq
            seq_lengths.append(len(seq))
            
            # 提取11到15列的数值
            mapping_quality.append(float(columns[10]))
            shift_dacs_to_pa.append(float(columns[11]))
            scale_dacs_to_pa.append(float(columns[12]))
            shift_pa_to_norm.append(float(columns[13]))
            scale_pa_to_norm.append(float(columns[14]))

    # 计算总体均值和标准差
    def calculate_mean_and_std(values):
        return np.mean(values), np.std(values)

    signal_mean, signal_std = calculate_mean_and_std(signal_means)
    single_base_mean, single_base_std = calculate_mean_and_std(single_base_signal_lengths)
    seq_mean, seq_std = calculate_mean_and_std(seq_lengths)
    mq_mean, mq_std = calculate_mean_and_std(mapping_quality)
    shift_dacs_mean, shift_dacs_std = calculate_mean_and_std(shift_dacs_to_pa)
    scale_dacs_mean, scale_dacs_std = calculate_mean_and_std(scale_dacs_to_pa)
    shift_pa_mean, shift_pa_std = calculate_mean_and_std(shift_pa_to_norm)
    scale_pa_mean, scale_pa_std = calculate_mean_and_std(scale_pa_to_norm)

    return {
        "Signal Values": (signal_mean, signal_std),
        "Base Signal Count": (single_base_mean, single_base_std),
        "Seq Length": (seq_mean, seq_std),
        "mapping quality": (mq_mean, mq_std),
        "shift_dacs_to_pa": (shift_dacs_mean, shift_dacs_std),
        "scale_dacs_to_pa": (scale_dacs_mean, scale_dacs_std),
        "shift_pa_to_norm": (shift_pa_mean, shift_pa_std),
        "scale_pa_to_norm": (scale_pa_mean, scale_pa_std),
    }

def main():
    parser = argparse.ArgumentParser(description="Analyze norm_signals_text, seq, and columns 11-15")
    parser.add_argument('--file', type=str, required=True, help="Input file path")
    parser.add_argument('--num_lines', type=int, default=-1, help="Number of lines to process (-1 for all)")
    args = parser.parse_args()
    
    results = process_file(args.file, args.num_lines)
    
    print("Overall Statistics:")
    print(f"  Signal Values Mean: {results['Signal Values'][0]:.4f}, Signal Std: {results['Signal Values'][1]:.4f}")
    print(f"  Base Signal Count - Mean: {results['Base Signal Count'][0]:.4f}, Std: {results['Base Signal Count'][1]:.4f}")
    print(f"  Seq Length - Mean: {results['Seq Length'][0]:.4f}, Std: {results['Seq Length'][1]:.4f}")
    print(f"  Mapping Quality - Mean: {results['mapping quality'][0]:.4f}, Std: {results['mapping quality'][1]:.4f}")
    print(f"  Shift DACS to PA - Mean: {results['shift_dacs_to_pa'][0]:.4f}, Std: {results['shift_dacs_to_pa'][1]:.4f}")
    print(f"  Scale DACS to PA - Mean: {results['scale_dacs_to_pa'][0]:.4f}, Std: {results['scale_dacs_to_pa'][1]:.4f}")
    print(f"  Shift PA to Norm - Mean: {results['shift_pa_to_norm'][0]:.4f}, Std: {results['shift_pa_to_norm'][1]:.4f}")
    print(f"  Scale PA to Norm - Mean: {results['scale_pa_to_norm'][0]:.4f}, Std: {results['scale_pa_to_norm'][1]:.4f}")

if __name__ == "__main__":
    main()
